#------------------------------------------------------------------------------
rm(list=ls())
library(LalRUtils)
libreq(tidyverse, data.table, zoo, tictoc, fst, 
       fixest, PanelMatch, patchwork,
       rio, magrittr, janitor, did, panelView, 
       ggiplot, tictoc, binsreg, interflex)
set.seed(42)
theme_set(lal_plot_theme())

#------------------------------------------------------------------------------



#------------------------------------------------------------------------------
# Define Paths
#------------------------------------------------------------------------------

# R studio
setwd( dirname(rstudioapi::getActiveDocumentContext()$path) )
# R default : unccoment if you use default R
# setwd(getSrcDirectory(function(){})[1])
#------------------------------------------------------------------------------


#------------------------------------------------------------------------------
# Load Data
#------------------------------------------------------------------------------

vcf <- fread("vcf_data_complete.csv", sep = ",")
gfc <- fread("gfc_dta.csv")
setnames( vcf, "d", "D")

#------------------------------------------------------------------------------


#------------------------------------------------------------------------------
# Panel A
#------------------------------------------------------------------------------

gfc = gfc[pref == 1] # subset to primary sample
# %%
xcreg = gfc[pesa_exposure == 1 & lag(pesa_exposure) == 0 ]
xcreg = xcreg[, dist2 := min_dist_to_mine * 100][dist2 <= 100]
# %% mining summary figure - binned scatterplot - fig8 panel A ;
gen = binsreg(y = xcreg$def_ha, 
              x = xcreg$dist2,
              nbins = 30, polyreg = 2)
f = gen$bins_plot + lal_plot_theme() +
  labs(y = 'Deforested Area in \nLast Year Before PESA', 
       x = 'Distance (km)') +
  theme( legend.position = 'none', 
         axis.title = element_text( size = 22 ), 
         axis.text = element_text(size = 20))

# %%
ggsave( 'appendix_figurea6PanelA.pdf', 
       f, device = cairo_pdf)

#------------------------------------------------------------------------------



#------------------------------------------------------------------------------
# Panel B
#------------------------------------------------------------------------------


# %% data prep for interflex
gfc[, `:=`(id_f   = as.factor(village),
           sty_f  = as.factor(styear),
           D_f    = as.factor(D)
)]
gfc[, distance := min_dist_to_mine * 110]
gfc[, mines_count_fac := as.factor(case_when(
  mines_in_5k %in% c(1,2) ~ "1-2",
  mines_in_5k %in% c(3,4) ~ "3-4",
  mines_in_5k > 5 ~ "5+",
  TRUE ~ "0"))]
#%% binning estimator
mining_all = interflex(Y = "def_ha", D = "D_f", X = "distance",
                       FE = c("id_f", "sty_f"), cl = "id_f",
                       data = gfc[distance <= 100 & pref == 1], theme.bw = T,
                       Xlabel = 'Distance', Dlabel = "PESA",
                       Ylabel = '\nForest Index (VCF)\n',
                       cutoffs = seq(0, 100, 10),
                       estimator = 'binning', CI = TRUE, 
                       cex.lab = 1.3,
                       cex.axis = 1.3, 
)
mining_all$figure
# %%
ggsave( "appendix_figurea6PanelB.pdf", mining_all$figure,
       device = cairo_pdf)

#------------------------------------------------------------------------------